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1 . 


INTRODUCTION 


This report consits of a technical explanation of the two 
computer codes that are employed to calculate two-dimensional 
inlet flow fields in a supersonic free stream. The first code 
(the Euler code) is a modified version of the program that was 
developed and used in reference 1 to calculate external flows 
about two-dimensional airfoils. The second is a modified version 
of a nonorthogonal grid-generation routine which was developed 
in reference 2 to generate body fitted coordinate systems for 
airfoils . 

In the next section of this report the Euler code is ex- 
plained in detail. Mathematical formulation of the problem is 
given. The computational algorithm and program logic are ex- 
plained. The third section of the report consists of a detailed 
explanation of the grid-generation routine. 


2. CALCULATION OF TWO-DIMENSIONAL INLET FLOW FIELDS 

(THE EULER CODE) 


In this section the technical approach employed for the 
calculation of two-dimensional inlet flow fields in a uniform 
supersonic free stream is presented. In Section 2.1 the mathe- 
matical formulation of the problem is given. Section 2.2 con- 
sists of the overall logic of the computer code that is used for 
the numerical integration procedure. Program usage and operation 
instructions are given in Section 2.3 and an explanation of pro- 
gram accuracy and limitations is given in Section 2.4. 



2.1. Formulation of the Problem 


In this section the governing equations (Euler equations) 
for inviscid flow fields in two-dimensional inlets are given, 
boundary and initial conditions are discussed and the numerical 
solution procedure and its stability are described. Our descrip- 
tion follows closely that of reference 1. 

2.1.1 Governing Equations .- The Euler equations in Carte- 
sian coordinates for two-dimensional or axisymmetric flows can be 
written in dimensionless, conservation-law form for a perfect gas 
without external forces as 2,4,3 



and 

p = (Y-l) fe - 0.5 p (u 2 + v 2 )] (2.2) 

Also, j = 0 for plane flows and j = 1 for axisymmetric flows. 

In these equations u and v are the velocities along the 
x and y coordinates respectively, p is the pressure, p is 
the density and e is the total energy per unit mass. The sound 
speed, a, is given by a 2 = y(y-l) - - 0.5 (u 2 + v 2 )J . In order 
to use a body-fitted coordinate system, the governing equations 
are re-written subject to the general transformation 
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£ = £ (x , y , t) 




n = n (x,y , t) > 
t = t 


(2.3) 


where E, is the coordinate along the body and n is the coor- 
dinate perpendicular to the body. 

In terms of the new independent variables E, , n and t 
which map the original x,y,t plane into a rectangular domain, the 
conservation form of the governing equations is 


8 q + 3 ~ E + 3 F + j H = 0 

t £ n J 


where 


= 2 

J' 


q = -, E = £ q + E E + £ r F / J , 


x 


f =(n t q + n x e + n y f )/j, h = - 


Here J is the transformation Jacobian 


J = £ H - E, n = 
s x y ^y x 


/ ( x r y^ - x y r 
/ 1 C i n £ 


(2.4) 


(2.5) 


The metrics £ ' , E, etc. are formed from the derivatives of 
x T , x^, etc., using the relations 


£ = Jy 

^x 2 n 


E, = -Jx 

y n 


n = -Jy_ 
x 2 E, 


n Jx. 
y = c, 


( 2 . 6 ) 


v E,, - -x E, - y C n. = -x ri - y ri . 

^ t t x 2 t y t t x 2 t y 
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We now define velocities along the E, and n coordinates 


U=£. + £u+£v 


(2.7) 


V + ^ + 11 u + n v 
t x y 

/% /\ 

Using these defined velocities, E and F can be written as 


/pU 


pV 


E = J 


- 1 


PuU + ^ x p 
pvU + £ y P 
' (e + p)U - 


,F = J 


-V 


puV + n p 
x 

pvv + 

(e + p) V - n fc p 


( 2 . 8 ) 


H = - 


/ PV 

-i / puv 
2 


PV 

(e + p) v 


The inlet flow solution is obtained as the steady state 
solution of the time-marching method employed to solve equation 
(2.4) . 


2.1.2 Boundary and Initial Conditions .- In the problem 
under consideration, four types of boundaries are considered 
(figure 1). These are solid boundaries, inflow boundaries, outer 
boundaries and outflow boundaries. 

Along the cowl surface and the ramp surface the tangency 
condition is satisfied, i.e. V = 0. This gives 
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( 2 . 9 ) 





The pressure on the body surface is found from a relation which 

is obtained by adding n times the ^-momentum equation to n. 

X 

times .the n -momentum equation. This gives 


d n + u9 n + v9 n 
t t t x x y 


- PU [ n x u c + n y v c 


_ 1EY 
y 


un + vn 
x y 


( 2 . 10 ) 


n £ + n K 

xx y y 


+ P, 


o 2 + n 2 

x y 


In this way the values of the flow field quantities on the wall 
are calculated by extrapolation from the inner field. 

Free stream values are specified at the inflow and outer 
boundaries. For economy in the use of mesh points in the work 
done to date, the inflow boundary was located a short distance 
downstream of the leading edge of the ramp. This requires special 
treatment for the mesh points on the inflow boundary which are 
downstream of the ramp-leading-edge shock. This treatment is 
discussed below under the paragraph dealing with initial condi- 
tions . 

The outflow boundary conditions on the portion of the out- 
flow boundary external to the inlet are calculated by zeroth- 
order extrapolation from the interior. On that portion inside 
the inlet duct, the boundary conditions are specified according 
to whether the condition being calculated is supercritical or 
subcritical. In the case of supercritical operation, the flow 
field variables are again calculated by zeroth-order extrapola- 
tion from the interior. For subcritical operation (subsonic 
outflow) a set of boundary conditions are incorporated which are 
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obtained from the steady-state forms of the governing equations 
assuming uniform parallel outflow. In accord with this assumption, 
a constant-area section is added to the downstream portion of the 
inlet duct and the axial velocity u, at the end of this duct ex- 
tension is fixed as constant at its value calculated from one- 
dimensional compressible flow equations. The velocity component 
parallel to the outflow boundary, v, is set to zero; pressure, p, 
and density, p, are calculated by zeroth-order extrapolation from 
the interior. Energy is calculated from these values of u, v, p, 
and p. For compatibility with the normal momentum equation, the 
pressure gradient in the direction perpendicular to the main flow 
direction, r) , is set to zero at both walls. 

The initial conditions are specified by using either impul- 
sive initial condition, i.e. freestream conditions or the final 
solution of a previously calculated flow field. When starting 
from freestream conditions, values calculated from 2-D shock 
theory are used at the 3-5 points on the inflow boundary that are 
downstream of the ramp shock. 

2.1.3 Solution Procedure and the Numerical Scheme .- In the 
solution procedure employed to integrate equation (2.4) a temporal 
linearization process is employed. In order to employ this on 
the vector terms of equation (2.4) one needs to evaluate the Jaco- 

A A A A A A A 

bian matrices A = 9E/9q and B = 9F/9q. The flux vectors E 
and F are both linear combinations of q, E, and F (reference 1) . 

In the solution procedure the linearized equations are cast 

6 7 

into delta-form algorithm by approximate factorization ' (AF). 

The algorithm is non-iterative, and requires the inversion of 
two block-tr idiagonal (4x4) coefficient matrices at each time 
step in the integration procedure at only the interior points. 
References 6, 7, 8 and 9 contain detailed descriptions and vari- 
ous applications of the delta-form AF algorithm, hence here we 
outline it briefly. The delta-form AF algorithm can be used 
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either with trapezoidal or Euler temporal implicit differencing 
and reads 


A A A , 1 A A A A 

(I + h6^A n ) (I + ) (q n - q n ) = - At(S^E + 6^F n + n H n ) 


( 2 . 11 ) 


-a At (6 rrrr + c5 ) q 
££££ nnnn ^ 


n 


Here 6_ and 5 are second-order central-difference operators, 

K n 

At is the integration step size, h = At or At/2 for first- 
order or second-order two-level time differencing, respectively. 
The metric coefficients are evaluated at time-level n + 1/2 for 
second-order accuracy. The fourth-order smoothing term is added 
to the RHS of the difference equation both to overcome non-linear 
instability and to damp short wave lengths (references 6 and 7) . 
It is generally desirable to add the smoothing term both in terms 
of explicit and implicit portions. The explicit part suppresses 
non-linear instabilities whereas the implicit part enables the 
use of a ~ O(At) for_ very large At, where a is the artificial 
viscosity coefficient. In the current work (reference 1) these 
terms have been included in a pseudo-implicit manner by inverting 
a product of scalar pentadiagonals 


-aAt(S 5?55 + « nnnn )q n = q* - q 11 (2-12) 

where 

[l + aAtJ 1 ( V^) 2 jJ + ctAtJ -1 ( V A ) 2 j] q* - q n = 0 

It should be noted that for compatibility with the boundary 
conditions, a second-order smoothing is used at points adjacent 
to the boundary, whereas no smoothing is used at the boundary 
points themselves. 
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2.1.4 Numerical Stability .- In explicit methods the step 
size in the marching direction is bounded by the Courant- 
Friedrics-Lewy (CFL) condition. According to this, for the AF 
scheme in the t - £ sweep 


At 


(O < i 

A max — 


(2.13a) 


Where (al) is the maximum spectral radius of the local 

a max 's 

eigenvalues of A. Similarly in the t-ri sweep, one could write 
the CFL condition as 


At 

An 


(0 Pmax - 1 


(2.13b) 


where ( o'! ) 

B max 

eigenvalues of 
problem can be 


is the maximum spectral radius of the local 

A 

B. A practical constraint in a two-dimensional 
written as 


At 


max of 



< 1 


(2.14) 


In implicit schemes usually this condition is relaxed so that 
the right-hand side constant can assume values much larger than 
one. In the present case typical values were around 5 to 10. 

But since the present scheme is only neutrally stable, its con- 
vergence cannot be improved by the use of very large time steps 
(reference 1) . 

It should be noted that in most of the calculations done to 
date (reference 13), a Courant number of about 10 was used with- 
out any stability problems. However the following precautionary 
remarks are in order: 
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(1) The cowl had to be introduced into the calculation over 
some 50 time steps by slowly inflating it (i.e. by slowly applying 
the wall boundary conditions) . 

(2) In the case of an inlet with an isentropic-compression 
ramp, calculations had to be made for about 500 time steps at a 
Courant number of one to capture the ramp shock alone. The cowl 
boundary condition was then introduced slowly as above and the 
Courant number was increased to 5 and then to 10. 

(3) In the case of subcritical operation, the calculations 
produced instabilities after the shock was expelled from the 
inlet. Stability was restored by keeping the Courant number 
equal to one. 


2.2 Overall Program Logic 

In this section the overall logic of the computer program 
used to solve equation (.2.11) is summarized. The details are 
given in reference 1. The computational procedure starts by 
forming the right-hand side terms (RHS) in equation (2.11), first 
the smoothing operator then the steady part. Then, the computa- 
tional algorithm forms Aq* by block tridiagonal matrix inversion 
in E, 


(I + h6^A n ) Aq* = RHS 


(2.15) 


Aq* is stored temporarily in q . The procedure continues 
with block tridiagonal matrix inversion in n 


(I + h6 B n ) Aq 11 = Aq* 

n 


Finally 


"n+ i 

q 


is obtained from 


(2.16) 
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i 



( 2 . 17 ) 


q 


n + 1 



A "n 

q 


Note that for the two-level scheme used here, the algorithm 
requires two levels of data storage. 

In the computer program used in this work the execution of 
the above outlined algorithm is performed by PROGRAM MAIN. It 
starts the computation by calling SUB.INITIA. This sub-routine 
defines the physical and mathematical constants that are to be 
used in the computation procedure. The initial field is also 
specified in this subprogram. An option is provided for starting 
capabilities from a previously calculated (converged) field. 

This subroutine also calls SUB. GRID by which either a simple- 
stretched grid for a flat plate is computed or a previously calcu- 
lated general grid is read from a disc. Then SUB.XYMETS is called 
to evaluate the metric coefficients for the grid system. 

The main program then calls SUB. EIGEN, which calculates the 
eigenvalues of the coefficient matrices A and B, it finds the 
spectral radius and finally computes the time-step from the pre- 
scribed Courant number. The optional SUB. MAP maps the P -matrix? 
this can either be the Jacobian of transformation or the calcu- 
lated pressure field. The normalized variables and P/P^ of the 
initial field, are then printed by SUB.OUTPT. Free stream quanti- 
ties are given or computed by SUB.HARVIO which is a library sub- 
routine . 

After the initial field is prescribed and printed, the main 
program starts to execute the integration loop. This is done by 
calling SUB. STEP. This subroutine in turn calls SUB.3C to calcu- 
late data at the flow field boundaries. SUB. STEP then calls SUB. 
RHS which computes the steady part of the difference algorithm 
and SUB. SMOOTH which forms the smoothing operator. The temporal 
differences are evaluated by SUB. DIFFER. SUB. STEP then calls 
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SUB.FILTRX which performs the block tridiagonal matrix inversion 
in the ^-direction and SUB.FILTRY which performs the block tri- 
diagonal matrix inversion in the rt-direction . Final output is 
furnished by subroutine SUB.HARVIO at the end of the iteration 
loop . 

It should be noted that the computer code employed here is 
an adaptation of a code which has been applied effectively to 
external flows over airfoils (reference 1) . The original code 
contains options for viscous effects (laminar and turbulent) , how- 
ever since no attempt has been made to include viscous effects in 
the present calculations, the pertaining subroutines have not been 
included into the above described overall flow logic. The origins 
program also contains options for upwind differencing in conserv- 
ative form (references 1, 5, and 6) to be used in transonic flow 
calculations. These subroutines have also not been used in the 
present work. 


2.3 Program Use and Operation 

In this section we give a detailed description of the input 
information required to run a typical case. As described in 
Section 2.2, the initial field as well as values of various con- 
stants are evaluated in SUB.INITIA. This subroutine also reads 
the input data except for the last card, which specifies geometry 
(planar or axisymmetric) and type of outflow boundary. This last 
card is read by SUB.BC. We now list the input data cards in 
order : 


Card No. Format 


Variables 


1 


815 


NMAX: 
JMAX: 
KMAX : 
NP : 


Maximum number of time steps 
Maximum number of points along E, 
Maximum number of points along n 
Number of time steps for calling 


SUB .MAP . 
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Card No . 


Format 


Variables 


1 


815 


METH: 


(continued) 


I RE AD: 


INVIS : 


2 


715 


IREGO: 


3 


4 


I STORE : 

JTAIL1 : 
JTAIL2 : 


315 


IPLOT : 
IUPWIND : 


IOSCIL: 


8F10. 0 


LAMIN : 

CNBR: 

DX: 


DY: 
FSMACH : 
SMU : 


Flag for upwind differencing. If 
METH > 0 , program skips upwind 
differencing . 

Option for grid system. If IREAD 
> 0, reads grid from disc. 

If INVIS > 0, program accounts for 
viscous effects; If INVIS = 0 , the 
calculation is inviscid. 

If IREGO > 0/ reads initial field 
from disc, otherwise starts from 
free stream conditions. 

If ISTORE > 0, solution data 
stored on disc . 

First ^-interior point, set to 1. 
Last ^-interior point, set to 
JMAX-1 . 

Index for calling plot routine. 

If IUPWIND > 0, skips upwind dif- 
ferencing . 

For stationary airfoil set IOSCIL 
= 0 . 

If LAMIN 0, calculates turbulent 
viscosity . 

Courant number; set to about 10. 
x-increment used in SUB. GRID. For 
grid read off the disc set DX = 0 . 
y-increment, as above. 

Free stream Mach number. 

Pseudo viscosity coefficient. Set 
to about 10 times the time step. 
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Card No. 


Format 


Variables 


4 

8F10. 0 

(continued) 

5 

4EI0.0 

6 

212 


EPS : 
RE: 
ALPHA: 

xoscilT 

VARA: 

VARB: 

VARC: 

LFAC : 


JAXI : 


Used in SUB. GRID. Here set EPS = 
Reynolds number 
Angle of attack 

Variables required for oscillating 
airfoil. For the inlet flow set 
to zero. 

If LFAC = 1 outflow is subsonic, 
otherwise it is supersonic 
If JAXI = 0 flow is plane 2-D, 

If JAXI = 1 flow is axisymmetric 


0 . 


In starting from free stream conditions, the number of points 
on the inflow boundary to be placed downstream of the ramp shock 
are prescribed in SUB.INITIA. Values for the variables are also 
prescribed in this subroutine using the 2-D wedge theory. 

The typical output of the program involves the following 
(described in more detail in reference 13): 

(1) Tabulation of the input parameters 

(2) Printout of the Jacobian matrix as calculated from the 
generated grid 

(3) A map of the transformation Jacobian 

(4) Printout of the flow field variables on the cowl surface 

(5) Free stream quantities given by SUB.HARVIO. Note that 
this is a library subroutine 

(6) Iteration index, physical time and time step 

(7) Maximum ERRORS occurring within the field 

(8) Values of LFAC (if LFAC > 1 outflow is subsonic) and 
JAXI (if JAXI = 0 flow is plane 2-D, if JAXI = 1 flow 
is axisymmetric) 
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(9) Mass flow rate evaluated at various ^-stations. 

(10) Residual quantity; must converge to a small number such 

— 1 — 2 

as 0(10 ) - 0(10 ) for converged solution 

(11) Line-printer plot of pressure in the computational plane 

(12) Line-printer plot of Cp on the ramp surface 

(13) Line-printer plot of Cp on the cowl surface 

(14) At time step equal to NMAX total output as furnished by 
HARVIO 


2.4 Program Accuracy and Limitations 

2.4.1 Accuracy . - The local accuracy of the finite difference 
scheme that is employed in this program is second-order both in 
time and space. It should be noted that the method is a conserva- 
tive difference formulation, therefore its global numerical accu- 
racy is also given by the local accuracy and is second-order in 
time and space. There are, however, other points that contribute 
to the accuracy of the numerical solution: 

(a) The grid system generated has a very important effect on 
the overall solution. With a course grid the shocks (discontinui- 
ties) are not adequately resolved, and in certain cases the re- 
sulting field is totally erroneous. Near the solid walls the 
grid system must be properly clustered. Grids also should vary 
smoothly in the flow field. 

(b) The overall accuracy also depends on the accuracy at 

the boundaries. In the present program at the flow field bound- 
aries it is assumed that q equals q . Therefore on the 

boundaries the formal accuracy is first-order in time. Second 
order accuracy in space, however, is retained all through the 
flow field. 

2.4.2 Limitations . - The main limitations of the program are 

(a) The flow field analyzed is two-dimensional. 


14 



(b) Owing to the neutral stability criteria, convergence to 
steady state cannot be accelerated by using large time- 
steps . 

(c) As currently implemented, the flow field analyzed is 
inviscid . 


3. GENERATION OF NONORTHOGONAL GRIDS FOR TWO-DIMENSIONAL INLET 
FLOW FIELDS 

In this section we present the method that is used to create 
nonorthogonal computational grids for the calculation of two- 
dimensional inlet flow fields. The computer program used is an 
adapted version of that described in detail in reference 2. There- 
fore here, we briefly present the outline of the mathematical 
formulation and the structure of the program. Detailed informa- 
tion is given only about the changes that were employed for the 
present purposes. 


3.1 Problem Formulation 

The computer program is based on the method developed in 
references 10 and 11 (TTM-method) to generate two-dimensional 
computational grids about arbitrary bodies with coordinate lines 
coincident with all boundaries. A provision is included for 
clustering about any grid line. 

The TTM method solves Laplace equations in the physical 

plane 


+ 

?yy 

= 0 

(3.1) 

n + 

'xx 

V 

= 0 

(3.2) 

to generate curvilinear coordinates 

£ (x,y) 

and g(x,y), which 
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form a rectangular grid in the transformed plane. In the compu- 
tational plane, equations (3.1) and (3.2) transform as 



ax £C 

- 23x. 

+ yx 

nn 

= 0 

(3.3) 


ay 5C 

’ 2Sy Cn 

+ yy 
, - / rin 

= 0 

(3.4) 

where 






a 

’ x n + *n' 6 

= x C x n + 

y £ y n and Y = x l + y l 

(3.5) 

The solution 

of equations 

(3.3) and 

(3.4) 

for x and y 

on the 


rectangular £ , ri computational grid yields the solution for £,n 
in the physical plane if the transformation Jacobian is everywhere 
nonzero. The most attractive aspect of the TTM method is the 
ability to arbitrarily locate boundary points by specifying their 
locations in the x,y plane. 

The program has the option of providing grid clustering along 
the y-direction to increase computational accuracy in regions 
where large gradients of the flow variables occur. This is espe- 
cially important in the present problem due to the existence of 
compression and expansion corners internal and external to the 
inlet. In this work, it was necessary to have the grid points 
clustered near the cowl surface and the ramp surface to insure 
accuracy. To this end the clustering transformation used in refer- 
ence 12 was employed; according to this transformation 


z 


a + ( 1 


g + z(2a + 1)- 2a " 
B - z ( 2a. + 1 ) + 2a. 


£n 


' 3 + 1 ' 
13 - 1 J 


(3.6) 
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permits the mesh to be clustered both near the cowl and the ramp 
surfaces when a = 1/2. Note that in equation (3.6) z is the 
spacing of constant £ lines in the computational plane, whereas 
z is the spacing of constant £ lines in the physical plane. 
Values of 3 are chosen for each constant-£ line from the re- 
quirement that z = z using a Newton-Raphson iterative 

max max c 

procedure. With 3 defined as the value of 3 after the nth 

n 

iteration 


3=3 
n n-i 


F(i W 


where 


F ( 3 ) 


max 


- z 


max 


(3.7) 


(3.8) 


Note that this procedure requires an initial guess for 3 Q on 
the first constant £- line. 

3.2 Overall Program Logic 

The execution of the foregoing calculation procedure is 
performed by PROGRAM MAIN. It starts the calculation by reading 
input data cards. SUB . INNER is called to distribute points on 
the rear boundaries and SUB. OUTER locates points on the bottom- 
front-top boundaries. Note that in our modified program all the 
bottom-front-top boundaries are read as input in SUB. OUTER. 

Initial conditions for the TTM are provided by equally spacing 
grid nodes along straight lines of constant £ between corre- 
sponding boundary points. The TTM method is applied by SUB. RELAX, 
stretching is applied by SUB.CLUSTR. Newton-Raphson scheme for 
finding 3 is applied by SUB.EPSIL, whereas SUB. TAINT, which is 
a system supplied subroutine, is used for polynomial interpolation 
of the tabulated values of the x^ , y^ coordinates as functions 
of z^. Note that the clustering procedure permits one to discard 
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the lines of constant n obtained from the TTM method and to 
recalculate them according to the desired clustering near the 
boundaries . 

Subroutines PLAWT, AXIS, LABEL, PLOT, and TITLE are used to 
illustrate the grids on the SC4020 plotter. They contain calls 
to the standard SC4020 software subroutines; thus, the SC4020 
software subroutines must be made available to the loader. If 
the user's computer installation does not include an SC4020, the 
five subroutines READIN, PLAWT, and EOF TV should be removed from 
the main program, and data card 7 should be deleted. 

3.3 Program Use and Operation 

In this section we give a description of the input information 
required to run a typical test case. Various subroutines have 
been modified to accomodate a flexibility for input data for the 
present calculations. We first list the input cards in order and 
then summarize the modifications. 

Input data is read by PROGRAM MAIN in the following order: 


Card No. Format 


Variables 


1 


315 


JWAKE : 


KMAX: 


MAX IT : 


Number of grid points in the 
x-direction in the wake includ- 
ing the trailing edge point. 

Set to 1 for the present calcu- 
lations . 

Number of grid points in the 
y-direct ion . 

Number of iterations to be used 
in creating the unclustered grid. 


2 


8F1 0 . 0 


XGMX : x-direction coordinate of rear- 
ward boundary. 

XGMN : x-direction coordinate of front 
boundary . 
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Card No. 


Format 


Variables 


2 

(continued) 


3 

4 


5 


6 

7 


8F10.0 


F10.0 
4F10 . 0 


F10 . 0 

F10.0 
A8 0 


YMAX: 
YMIN : 
XNOSE : 
XTAIL s 
XRAMP : 

DY1 : 

XORG: 
YORG : 
ETAC: 
BETA: 

OMEGA: 

DY2 : 
A8 0 : 


y-direction coordinate of front 
boundary. 

y-direction coordinate of bottom 
boundary. 

x-direction coordinate of leading 
edge of cowl. 

x-direction coordinate of end 
of cowl. 

x-direction coordinate of lead- 
ing edge of ramp . 

Minimum y-increment .on rearward 
boundary for initial conditions. 

Location on the x-axis of the 
special clustering origin. 
Location on the y-axis of the 
special clustering origin. 

Angle (in degrees) , about which 
angular clustering is done. 
Parameter determining the 
strength of angular distribution 
set to zero. Otherwise values 
in the range 1 to 5 . 

Parameter for Thompson solver 
values between 0 to 2 recom- 
mended . 

Minimum y-direction spacing for 
final clustering. 

Computer card for the title 
frame of SC4020 plots. To be 
deleted if the user’s 
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Card No. 


Format 


Variables 


7 


A8 0 


{continued 


8 


A71 


9 


A8 0 


10 


2F10.0 


installation does not include an 
SC4020 . 

Title on each SC4020 plots re- 
resulting from each run. 

Description of the cowl. 

All cards following card 9 con- 
tain in columns 1 to 10 and 11 
to 20 two floating-point numbers, 
x and y, respectively, which are 
coordinates of points on the cowl 
surface defining the cowl shape. 
The points should proceed counter- 
clockwise . 


It should be noted that the number of points in the £- 
direction is specified in PROGRAM. MAIN by assigning an integer 
value to the variable NBOD . We have further modified SUB. OUTER 
such that it allows the boundary points along the ramp inflow 
boundary and the outer boundary to be specified point by point. 
For each different geometry, these boundary points must be pre- 
scribed by making the required changes in SUB. OUTER. 

This completes of discussions on the overall aspects of the 
grid generation program that was employed to generate non-orthog- 
onal grids for two-dimensional inlet flow fields. Example calcu- 
lations to illustrate the use of the program are reported sepa- 
rately (reference 13) . 
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